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Abstract 

The ground state persistent current and electron addition spectrum in two- 
dimensional quantum dot arrays and one-dimensional quantum dot rings, 
pierced by an external magnetic flux, are investigated using the extended 
Hubbard model. The collective multidot problem is shown to map exactly 
into the strong field noninteracting finite-size Hofstadter butterfly problem 
at the spin polarization transition. The finite size Hofstadter problem is dis- 
cussed, and an analytical solution for limiting values of flux is obtained. In 
weak fields we predict novel flux periodic oscillations in the spin component 
along the quantization axis with a periodicity given by v h/e (y < 1). The 
sensitivity of the calculated persistent current to interaction and disorder is 
shown to reflect the intricacies of various Mott-Hubbard quantum phase tran- 
sitions in two- dimensional systems: the persistent current is suppressed in the 
antiferromagnetic Mott-insulating phase governed by intradot Coulomb inter- 
actions; the persistent current is maximized at the spin density wave - charge 
density wave transition driven by the nearest neighbor interdot interaction; 
the Mott-insulating phase persistent current is enhanced by the long-range 
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interdot interactions to its noninteracting value; the strong suppression of the 
noninteracting current in the presence of random disorder is seen only at large 
disorder strengths; at half filling even a relatively weak intradot Coulomb in- 
teraction enhances the disordered noninteracting system persistent current; in 
general, the suppression of the persistent current by disorder is less significant 

in the presence of the long-range interdot Coulomb interaction. 
PACS numbers: 73.20.Dx, 73.23.Hk, 71.10.Fd 
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I. INTRODUCTION 



In this paper we consider an array of coherently couple semiconductor quantum dots 
arranged in finite two-dimensional (2D) square lattices or one- dimensional (ID) rings pierced 
by a magnetic flux oriented normal to the lattice or the ring plane. At low temperatures, 
these quantum dot arrays may be considered "artificial molecules" (with individual quan- 
tum dots being the "atomic" constituents of these artificial molecules) because the electron 
phase coherence length is comparable to the array linear size. Theoretical work on mul- 
tidot systems has mostly concentrated on the two limiting situations: coherent dots with 
no Coulomb interaction! and Coulomb blockade of individual dotsii. In this paper, we 
consider quantum dot arrays taking into account quantum fluctuations arising from inter- 
dot hopping, electron-electron interaction, and random disorder effects through an extended 
Hubbard-type HamiltonianUl^ 

In a previous paper ,0 we reported on our prediction of an equilibrium persistent current 
in finite 2D dot arrays (without any periodic boundary conditions) in the presence of an 
applied magnetic field transverse to the 2D plane. In this paper we provide details and 
expand on our previous work, and present results for the electron addition spectrum and 
the persistent current in 2D square lattices and ID rings including effects of collective physics 
arising from the multiple dot structure of the system within a simple model for the single- 



understand lattice effects on the persistent current and the electron addition spectrum, the 
lattice here being the artificial lattice defining the 2D or the ID quantum dot array (with 
the typical lattice constant in the 20 — 200 nm range). One goal is to identify experimentally 
observable features of various Mott-Hubbard quantum phase transitions (including realistic 
disorder and interaction effects) in semiconductor quantum dot arrays. 

The importance of electron-electron interaction has been stressedEUll in the literature 
in the context of persistent current experiments in ID gold and semiconductor ring si. The 
magnitude of the persistent current^ in disordered gold rings was found to be one to two 




One of our primary motivations is to 
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orders of magnitude larger than that theoretically predicted, whereas in clean semiconductor 
rings the magnitude of the persistent current was found to be in a good agreement with the 
theoretically predicted simple noninteracting value of e vp/L (with vp being the Fermi 
velocity of electrons moving in a ring of length L, in our notation L denotes the size of 
the system as defined by the total number of dots in it). Although the interplay between 
disorder and Coulomb interaction in determining the magnitude of the persistent current 
in ring topologies is the subject of many recent theoretical investig ationsElffl, the issue 
remains unsettled. Disorder and interaction effects are naturally included in our Mott- 
Hubbard model of finite quantum dot arrays, and we will comment on their influence on the 
persistent current. 

The paper is organized as follows. In Sec. II we define and describe the extended Mott- 
Hubbard Hamiltoniani which forms the basis of our theoretical description of the collective 
physics in finite quantum dot lattices. In Sec. Ill we present our calculated electron addition 
spectrum as a function of the externally applied magnetic flux for a 3 x 3 lattice and an 
L = 9— site ring, and identify the main features of these results which are studied in the 
subsequent sections. In Sec. IV we clarify the physical meaning of the magnetic field depen- 
dence of the addition spectrum by demonstrating the equivalence between the derivative of 
the Hamiltonian with respect to the flux and the magnetization density (or equivalently the 
persistent current) operators. In Sec. V we study the energy spectrum and the persistent 
current of the finite open 2D \fh x \fL tight-binding lattices of L noninteracting quantum 
dots. We also discuss in this context the Hofstadter spectrum of an infinite tight-binding 
lattice. We establish a connection between the ground state magnetization in the lattice 
and in the continuum 2D system by identifying the different regions in the energy spectrum 
in a lattice, and illustrate them with the calculated distributions of the persistent current 
on a 15 x 15 lattice. We also solve the problem exactly for lattices of arbitrary sizes in 
the two limiting situations of the flux = and = 0.5 0o per unit cell (with 0o = h/ e 
being the fundamental flux unit), and classify the = states using group theory and 
perturbation theory. We also discuss in Sec. V our results for the noninteracting persistent 
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current in ID rings. In Sec. VI we study the electron addition spectra of interacting lattices 
in a magnetic field. In Sec. VIA we study intradot Coulomb interaction effects within the 
minimal Hubbard model approximation. In particular, we study the periodic oscillations of 
the component S z of the total ground state electron spin along the quantization (z) axis of 
the 2D Hubbard model using the Lanczos exact diagonalization technique, and by solving 
the Bethe ansatz equations in ID rings. We reanalyze the persistent current results in a ID 
Hubbard ring, and find that the interacting system behaves as a single particle by changing 
its total orbital momentum sequentially as a function of the flux. The details of our Bethe 
ansatz analysis for the ID Hubbard ring spectrum are given in Appendix A. In Sec. VIA 
we also discuss the finite size realization of the Mott-Hubbard metal-insulator transition 
and magnetic ordering in finite 2x2, 3x2, 4x2 and 3x3 quantum dot clusters (with 
L = 4, 6, 8, 9 dots respectively in the system). In Sec. VIB we discuss the spin density 
wave-charge density wave ordering transition in a half-filled 2D 3x3 array in the presence of 
the nearest neighbor interactions. We find an enhancement of the persistent current along 
the transition line. In Sec. VII, we consider random disorder effects on the 2D persistent 
current. We find that the noninteracting persistent current as a function of the disorder 
strength shows a behavior similar to that of the conduct ivityS it is strongly suppressed 
only at large disorder strengths. We obtain an empirical scaling function for the persistent 
current as a function of the disorder strength. Finally, we discuss the effects of having both 
disorder and Coulomb interaction on the 2D persistent current in the 3x3 array by using 
the Lanczos exact diagonalization technique. We conclude with a summary of our results in 
Sec. VIII. 

II. MODEL 

We model an isolated finite system ("array") of coherently coupled (nominally identical) 
semiconductor quantum dots arranged in one-dimensional rings or two-dimensional rectan- 
gular ( "square" ) lattices at zero temperature. We assume that the charging of an otherwise 
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electrically neutral quantum dot array with a fixed number N of excess quasiparticles is 
accomplished through the tunneling of the quasiparticles from a nearby backgate electrode. 
The excess charges are shared among all dots in the array in a molecular-like fashion by 
quantum-mechanical tunneling, i.e., our quantum dot array is coherent. The equilibrium 
properties of quantum dot arrays we investigate in this work can in principle be measured 
through experiments on semiconductor dot arrays using available experimental techniques: 
tunneling transport spectroscopy^, equilibrium capacitance spectroscopy^, and equilibrium 
magnetization measurements^. 

The Hamiltonian of an isolated array of L coupled quantum dots is given by the sum of 
three terms: 

H array = (KE) sp + (KE) hop + V int . (1) 

The (KE) sp term in Eq. [I] includes all intradot single-particle effects (including con- 
finement contribution) and is the sum of the single-particle Hamiltonians of the individual 
quantum dots of the array. We model GaAs quantum dots at zero temperature. The char- 
acteristic size of each dot is taken to be D > where ag(~ 100 A) is the effective Bohr 
radius of the bulk GaAs material. All effects of the electron-electron interactions within 
a neutral dot are absorbed in the effective values of quasiparticle parameters giving the 
effective mass m* = 0.067 m e and the effective g-factor g* = 0.2 g e in a standard band 
description for the single-particle energy levels within the dot. The confinement potential 
of a quantum dot is known to be approximately parabolic^. Due to quantum confinement 
in the dot, a continuous conduction band for the excess quasiparticles is a discrete series of 
single-particle energy levels e a where a denotes a single-particle state including spin. The 
single-particle intradot level spacing at zero field is taken to be A = Hujq (with ujq essentially 
being a harmonic oscillator frequency). In our work we consider only the lowest intradot 
states near the Fermi level. In the occupation basis of states |a) the (KE) sp term is written 
in the second quantized notation as 
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L 

{KE) sp = £ia(B)c ] ia c ia . (2) 

In Eq. the summation is over all dots i\ c^i a (ci a ) is a creation (annihilation) operator for 
a quasiparticle on the zth dot in a state a. The single-particle magnetic field dependence is 
included in Eq. ^ through the usual Fock-Darwin-Zeeman scheme as 

e ia = e a = h[(uj c /2) 2 + u 2 ] 1/2 + {-l) a g^ B B/2. (3) 

with u c = eB/m* being the cyclotron frequency. We consider a single spin-split level per 
dot, and set a = 1, 2 in Eq. |3| to correspond to the spin up/down lowest confined quantum 
dot level. We neglect correlations arising from single-particle level crossings that have to be 
taken into account for fields larger than B = [g*(g* + 2m /m*)] -1 / 2 (A / jis) , and concentrate 
on the collective physics in the array. 

The quasiparticles are allowed to tunnel ("hop") between the same single-particle states 
| a) in the dots with the tunneling amplitudes t a = t. In the tight-binding approximation, 
keeping only nearest neighbor tunneling, the tunneling energy is given by 

{KE) hop = £ (t^cUaCja + h.C.), (4) 
<i,j>,a 

where 0^- = | ■ A ■ Uj is the Peierls phase factor^!, with A as the magnetic vector potential. 
The indices i, j denote the spatial positions of the dots. The (KE)h op term defines the 
topology of the array. We model a one-dimensional ring of L quantum dots with a total 
magnetic flux piercing its enclosed area and a two-dimensional square lattice of L = L x x L y 
with open boundary conditions with a magnetic flux piercing each unit cell of the lattice. 

The third term in Eq. (|I|) defines the intra- and interdot Coulomb interactions between 
quasiparticles 

ij 

where pi = J2a c ^ia c ia is the number operator for the zth quantum dot. The interaction 
constants are related to the capacitance matrix of the quantum dot array by Vij = (C" 1 )^-, 
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where Cu = C g + NiC and CV,- = —C for nearest neighbor dots.mOti The capacitance C g 
represents the capacitance of a quantum dot with respect to the external gates, while C 
denotes the capacitive coupling between the zth dot and the iV, neighboring dots. Eq. ([]]) 
thus has the form of an extended Hubbard model with screened long-range interactions. For 
C <C C g , the interaction matrix elements fall off as Vy ~ {7(C/C £ ,)' l ~ J 'L We include the 
effects of short range interactions, keeping only the on-site interaction U and the nearest 
neighbor interaction V. We use t = 0.1 meV, A = 3t, and U = lOt in our calculations (unless 
otherwise stated) as a representative set of Hubbard parameters describing the GaAs dot 
arrays. The Hamiltonian given by Eqs. [3] - |^ has been used earlier to describe coherent ID 
quantum dot chains with open or periodic boundary conditionsi'0. (Inclusion of disorder 
in our model will be discussed in Sec. VII.) We calculate the electron addition spectrum by 
doing an exact diagonalization of Eq. [I] in the subspace of the total number of quasiparticles 
N in the array and the total spin component S x (= —\{N — M) + \M with M being a number 
of spin-up electrons) along the external magnetic field B. The Hilbert space of Eq. |I| with 
fixed N and S z grows exponentially with the system size. We use the Lanczos methodii 
for our exact diagonalization of Eq. [TJ, and carry out a ground state energy minimization 
over S z to find the stable ground state for a given N. The largest matrix size that we have 
considered is 15876 2 which corresponds to the interacting L = 9 (or 3x3) dot array at 
half-filling for N = L = 9. We also calculate the energy spectrum of one-dimensional rings 
in the limit C = (Va = U ^ 0, V^- = 0) by solving numerically the corresponding Bethe 
ansatz exact solution equations. 

III. ELECTRON ADDITION SPECTRUM 

By definition the chemical potential /ijy of the array is given by 

m = E (N)-E (N-l), (6) 

where Eq(N) is the minimum eigenvalue (i.e. the ground state energy) of Eq. [I] in a space 
of fixed N and all allowed M. The addition spectrum of the system is the /ijv — N plot, 
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which we show in Figs. 1 and 2. The calculated chemical potential of a 3 x 3 array and a 
9— site ring in the minimal Hubbard model approximation (Vu — U ^ 0, Vij — 0) is shown 
as a function of the applied flux <p in Fig. 1 . Each curve in Fig. la (c) traces the chemical 
potential /z a? of a 3 x 3 array (L = 9 sites ring) with N electrons as a function of the magnetic 
flux (j)/(f)o through a unit cell (the ring), with O = h/e. The z-component (S z ) of the total 
spin of the corresponding ground state of the system is shown in Fig. 1 (b) and (d), and 
the critical magnetic flux for full spin polarization in the array and the ring is given as an 
inset in Fig. 1 (a) and (c), respectively. 

The three main features of the results shown in Fig. 1 are: (a) the chemical potential 
spectrum evolves with the maximization of the total spin polarization in the system; (b) 
apart from an aperiodic single-particle background contribution, fi(N) is a periodic function 
of the magnetic flux with a flux periodicity of 7 0o with 7 < 1; (c) the spin polarization 
transition of the system occurs through the cycles of the periodic S z oscillations in the weak 
magnetic field region. 

The result (a) is a trivial outcome of the minimal Hubbard model: a dominant Zee- 
man energy term g*fj,B leads to single-occupancy, effectively suppressing the presence of 
the Hubbard U term in Eq. |I|. The chemical potential spectrum therefore behaves as that 
of noninteracting spinless fermions in the maximum spin-polarization region. The energy 
spectra of the 9— site systems plotted in Fig. 2(a) and in Fig. 7 for noninteracting quasi- 
particles [setting = A, N = M in Eq. |I] can be directly compared with the strong field 
regions in Figs. 1(a) and (c). In the weak field region the single-particle physics can again 
be distinguished from the collective physics, and the spectrum in this regime is plotted in 
Figs. 2(b) and (c) for the array and the ring, respectively. In Figs. 2(b) and (c) we set 
£j = A, £-[■ = A + 5, where 5/A(= 0.03) <C 1, and M corresponds to the minimum eigenvalue 
of Eq. [I] for a given N. (We keep a small field-independent shift 5 between opposite spin 
single-particle states for our later discussion of the S z oscillations.) The parameters of the 
Hamiltonian that we use in Eqs. ([!])- (|§) determine the value of the critical magnetic flux <p c 
needed to produce full spin polarization in the 3x3 array to be <p c m 8t 2 a 2 / \g* [ibU) ~ 520 o 
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with a lattice constant of a = 280 nm. For illustrative purposes we rescaled the magnetic 
flux by 1/32 to show the full behavior of the spectrum on a single scale in Figs. 1(a) - (d). 
(This rescaling physically corresponds to rescaling the lattice constant to ~ a/6.) 



IV. PERSISTENT CURRENT 

We first clarify the physical meaning of the intricate magnetic field dependence of the 
energy spectrum shown in Figs. 1 and 2. It is well-known@ that gauge invariance along 
with the single-valuedness of the electron wavefunction allows for the existence of a ground 
state persistent current in normal metal rings threading an external magnetic flux. The 
intrinsic magnetic moment associated with this persistent current, which is proportional to 
the persistent current itself in ID rings, is an oscillatory function of the external flux with 
a period equal to the elementary flux quantum <p . The existence of such an oscillatory 
persistent current in normal metal rings has been experimentally verified!!. 

In finite 2D continuous systems no magnetization is expected classically.0 It is argued 
that the magnetization due to electron orbits along the edge of the sample exactly cancels the 
magnetization arising from the bulk orbits. In a quantum-mechanical description, however, 
the contribution from the edge states is expected to be statistically insignificantil, and the 
bulk contributions lead to the famous Landau diamagnetism in macroscopic 2D systems. In 
mesoscopic systems the phase coherence length is comparable to the linear system size 
L, and as was shown in several theoretical papersil'0, the edge states in this situation can 
carry a persistent current creating a paramagnetic moment in the continuous 2D geometries 
(e.g., a 2D disk-shaped quantum dot) comparable in magnitude to the Landau diamagnetic 
term. The Aharonov-Bohm effect leads to a flux periodicity of the magnetization carried by 
the edge states effectively forming a ID ring geometry in continuous 2D mesoscopic systems. 
We characterize the ground state magnetization of finite 2D quantum dot arrays without 
resorting to an artificial separation of bulk and edge states (which are not really meaningfully 
distinguishable in small structures) by considering the lattice model. The lattice spectrum 
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contains the "edge" states, the "bulk" states, and all other electron states given by the 
superposition of all topologically closed electron paths in the finite 2D lattice. 

The part of the Hamiltonian in Eq. [I] leading to the persistent current is the kinetic 
energy "hopping" (KE)h op term in Eq. |j. The total current operator^ is given by the 
commutator of the Hamiltonian with the polarization operator 

J T = ~ [H, P] , where P = -e ^ R^. (7) 

In a discrete lattice subjected to a magnetic field, the above formula becomes 

J t = ~ E ( R j - R i) Im K> ^ch a c ja } . (8) 

00 <ij>f a 

The expectation value of the total current operator is zero in the ground state of an isolated 
finite quantum dot lattice. Then each divergenceless term in the sum of Eq. |8] can be 
identified as the current between two nearest neighbor lattice sites i and j: 

Jjj — ^ y (-^j R«) Itti i t ij a e ^ ^ i a Cj a \ . (9) 

00 a k J 

The magnetic moment operator is 

M = - / R x J(R) d 3 R. (10) 



2 

The magnetization density operator has a non-vanishing z-component (i.e. along the mag- 
netic field direction) given by 

M. 1 



m 7 



Area 2 n cMs <{j> 



E J ij \. x i(yj - Vi) - yi( x j - x i)] ■ ( n ) 



Apart from the field dependence of single-particle levels in a single dot, the derivative of 
the Hamiltonian in Eq. p] with respect to the flux is 

OH dH 2 



E \^r) Im {%«^4c4. (12) 



<9$ n ceUs d0 n cMs . (! 

In a symmetric gauge, A = (— -B|, JB|, 0) in a uniform field, and one gets for the Aharonov- 
Bohm phase 0if 
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(13) 



The two operators in Eqs. |TT] and [12| are equivalent, and therefore 

dH 



m. 



(14) 



n C eiisd(p 

Eq. [14] is the usual thermodynamic expression for the magnetization, which we have derived 
here for our microscopic model. 

In this paper we use the convention of calling the magnetization density m z the persistent 
current I in both ID and 2D systems: 

/ = m z = (15) 

n C ellsO(P 

Experimentally, the equilibrium persistent current is usually observed by measuring the 
ground state magnetization. 

V. THE NON-INTERACTING SPECTRA 

We consider first the Hofstadter probleml! of a single particle in a magnetic field on an 
infinite tight-binding lattice. In the Landau gauge, A = (0, Bx, 0), the discrete Schrodinger 
equation in the occupation basis \a) = ^^(x, y) is 

fay) 

ip{x + 1, y) + ip{x — 1, y) + e t27T ^o x ip(x, y + 1) + e 2n ^o x ip(x, y — 1) = —E/t if)(x, y), (16) 

where <fi is the flux through a unit cell, and the lattice constant a is taken as the unit length. 
At = 0, an infinite system described by Eq. [16] is translationally invariant. The coefficients 
of Eq. [0] involve only x. The y— motion separates outS assuming that the y— part of the 
wavefunction preserves its <fi = plane wave form: 

1>(x,y)=A{k y )il>(k y ,x)e*y. (17) 

The function ifj(k y ,x) is a solution of 
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ip(x + 1) + ip(x - 1) + 2 cos 2tt^-x - fc y ^(ac) = -E/t tp(x). (18) 

Eq. |18| is the well-known Harper equation, describing a particle in a one-dimensional quasi- 
periodic potential. An infinite tight-binding lattice in a nonzero flux is no longer invariant 
under lattice translations. It has, however, been shown that it is invariant under magnetic 
translations^ by qa, whenever the flux through a unit cell is a rational fraction <f> = 2 (f> of 
the fundamental flux quantum with p and q being any two integers. The Harper equation 
spectrum has q energy bands at the rational fractional values of flux = 2 </>q. The magnetic 
translations by qa define a magnetic unit cell with the total flux = p<p through the 
magnetic unit cell. At incommensurate flux values (i.e. when 0/0o is n °t rational) the 
spectrum becomes a Cantor set.H The Cantor set spectrum has an infinite number of energy 
bands that exhibit a self-similar multifractal behavior. The energy spectrum of Eq. [18] is 
always a continuous function of the magnetic flux,il independent of whether 0/0o is rational 
or irrational. This enables one to make a direct connection between the lattice and the 
continuum spectra, for example, in identifying Landau bands in the lattice spectrum. 

The basic features of the spectrum!! of Eq. [18] (which we will refer to as the Hofstadter 
spectrum), which are also pertinent for a finite lattice, are: (a) it is periodic in <p Q , E((p) = 
E((p + n<fto); (b) it is an even function of flux, E(<p) = E(—<f))] (c) both E{(f>) and —E(<f)) 
belong to the spectrum; (d) the spectrum is bounded, —At < E(<ft) < At. A quarter of the 
Hofstadter spectrum plotted for rational values of the flux is shown in Fig. 3 where we 
explicitly label the first two Landau bands. The Landau bands are formed in the continuum 



limit when the magnetic length far exceeds the lattice constant, i.e. Iq = y (l/27r)(0 o /0) ^ 
1. In the tight-binding model, the effective mass is m ~ H 2 /2t. The energy levels in a lattice 
for l ^> 1 are approximately given by the continuous system Landau level expression, 
E n = hu c {n + |) with Tiuj c = Ant 0/0 o 0- For example, the ratios of the slopes of the first 
three Landau bands in Fig. 3 obey 23.4 : 15.1 : 5.5 w 5 : 3 : 1. The analytical solution of 
Eq. [18| for rational flux values were recently obtained using the Bethe ansatz method. El The 
complete spectrum is extremely complex, but a general feature of the spectrum which can 
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be seen in Fig. 3 and to which we will later return in our study of the persistent current in 
finite systems, is the presence of energy bands separated by large gaps. 

A quarter of the spectrum for a finite \[L x \[L = 15 x 15 (L = 225) lattice is shown 
in Fig. 4. A qualitative similarity between the spectra plotted in Figures 3 and 4 was 
pointed out in the literature^ the presence of similar energy bands in the spectrum, where 
the gaps between the bands are filled by the edge states which necessarily exist in finite 
systems. A finite spectrum was studied earlier in connection with the Quantum Hall effect 
and mesoscopic Aharonov-Bohm fluctuations. The emphasis of these earlier studies was on 
the part of the spectrum where both the Landau bands and the edge states can be clearly 
identified (the region from (d) to (e) in Fig. 4). Using Eqs. ([|) and @, we calculate the 
persistent current distributions in finite lattices for all eigenstates of Eq. [L6|. We follow Ref. 
|31| in the identification of different regions of the spectrum and illustrate each region with 
a sample current and a charge density distribution shown in Fig. 5. The lower-half band 
states in the weak field^ regime (Z > y/L) labeled (a), (b) and (c) in Figures 4 and 5 can 
be considered to be extended bulk states because the paths of the persistent current carried 
by these states extend across the sample with finite weights both at the boundaries and 
in the bulk. The sign of the current carried by these states oscillates, but as <fi — > it is 
determined by the degree of the degeneracy g of the spectrum at <fi = 0. We identify the first 
three Landau bands in the spectrum labeled (d), (e), and (f) in Figures 4 and 5 (the traces 
of the 4th and 5th bands can also be seen in the plot). The ratios of the slopes of the first 
three bands are comparable to those of the infinite tight-binding lattice, 20.2 : 14.1 : 5.3. It 
is interesting to note that the radii of the "large weight" persistent current orbits in Fig. 5 
approximately satisfy the semiclassical expression for the cyclotron radius of the nth Landau 
band, R n = lo\/2n + 1. The ratios of the radii of the persistent current orbits of the first 
three Landau bands in Fig. 5 obey 5:4:2^ \/E : \f?> : 1. Thus the continuous system 
result approximately holds for a finite lattice as well: each nth Landau orbit accommodates 
one more flux quantum than the (n — l)st band orbit. 
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For larger flux > f 2 -, the Landau levels form a more complicated but less degenerate 
pattern. We label the representative states (g), (h), and (i) for this region. The distribution 
of the persistent current carried by these states appears to consist of disconnected orbits, 
which can be found anywhere in the lattice. Intuitively, the existence of these disconnected 
persistent current orbits is consistent with our general expectation that smaller quantization 
orbits should not be strongly affected by the confining potential. 

The gaps of the infinite system in Fig. 3 are filled with edge states between the bulk 
Landau levels and the branched Landau levels in Fig. 4. The current distributions of 
the representative edge states labeled (j), (k), and (1) in Fig. 4 are shown in Fig. 5. 
These edge states have the largest weights concentrated near the boundaries in agreement 
with these states being called the "edge" states. The edge states in the regions between 
branched Landau levels (state (1) in Figures 4 and 5) seem to have a more complicated 
spatial distribution than the state labeled (j). The persistent current orbits of the state (1) 
are still extended across the lattice. The sign of the current carried by edge states of the 
type (j) is paramagnetic in accordance with the semiclassical argument given by Peierls@. 
The flux and energy separation between these states correspond to one added flux quantum 
through the total area of the sample. This is the origin of the periodic Aharonov-Bohm 
oscillations superimposed on a regular pattern of de Haas-van Alphen oscillations in the 
magnetization and in the magnetoconductance discussed earlier by Sivan and Imrylll. 

The total persistent current J of a system of N noninteracting spinless particles for two 
values of N selected in Figures 4 and 5 is shown in Fig. 6. From the symmetries of the 
L = L 1 / 2 x L 1 / 2 lattice spectrum, it follows that I{—<j>) = —I(4>), I(4> + n<t>o) = !{4>)i an d 
I(N) = I(L — N). The total current is given by the sum of currents carried by all occupied 
single-particle states, and it can be quite different from the persistent current carried at the 
Fermi energy. For example, in Fig. 6(a), at flux <p = 0.25 <po, the state at the Fermi energy is 
a 2nd Landau level bulk state (denoted (e) in Figures 4 and 5) which carries a diamagnetic 
persistent current. But the total persistent current is positive at this value of flux in Fig. 
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6(a), corresponding to the positive contributions from each filled edge state. The change 
in the sign of the current from paramagnetic to diamagnetic in the region from pa 0.1 0o 
to ~ 0.2 0o corresponds to merging all first 26 levels into the lowest Landau band. This 



merging of edge states into bulk states was discussed in Ref. |31|. A similar change of the 
sign of the persistent current can be seen in Fig. 6(b) for the region of flux from pa 0.3 0o 
to ~ 0.37 0o, where the Fermi Energy crosses from the branched Landau states through 
an edge state into the lower Landau level branched states. 

As mentioned above, an analytical solution for the eigenvalues of the Hofstadter problem 
(Eq. [18]) at rational values of flux on an infinite tight-binding lattice has very recently been 
obtained using the Bethe ansatz technique.il The fact that an infinite 2D problem (Eq. [16]) 
can be reduced to an effectively ID problem (Eq. [0|) is instrumental in using the Bethe 
ansatz technique. In a finite 2D lattice with open boundary conditions, a general solution 
is formally written as 

il){x,y) = A(k y ) ip(k y ,x) exp tkyV +A(-k y ) ip(-k y ,x) exp~ lkyV , (19) 

where ip(k y , x) and tp(—k y , x) are the formal solutions of Eq. [L8|. The motion in the x and y 
spatial directions generally can not be separated in Eq. (fL9D . The difficulty encountered in 
obtaining an analytic solution for the spectrum of a finite lattice stems from our inability to 
factor out a ip(x) factor in Eq. [19] corresponding to the opposite k y momenta of an infinite 
lattice. The origin of this is the breaking of the time-reversal symmetry by the magnetic 
field. Eq. [L8| is not invariant under the k y — > —k y transformation. However, for the two 
limiting values of the magnetic flux = and = 0.5 0o, the time-reversal symmetry is 
not broken, and the finite lattice problem can be analytically solved by reducing it to two 
effectively one-dimensional problems. The exact solution for these two limiting values of 
flux is given by 

ip(x, y) = A(k y ) ip(k y , x) sm(ky), (20) 

where k y takes on discrete values k y = (n = 1,2, ...,L y ) with L y being the extension 
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of the sample in the ^—direction. The x part of the wavef unction in Eq. (|20"D satisfies the 
Harper equation, 

ip(x + 1) + tp(x - 1) + 2 cos (^TT^-x^j cos(k y )i/j(x) = -E/t ip(x). (21) 

It can easily be verified by direct substitution that Eq. fl2HD is indeed a solution of Eq. ( |TED 
for the two limiting values of the flux, = and = 0.5 0o- 

At — 0, Eqs. ( pOl) and fl2l|) reduce to a trivially diagonalizable problem with a solution 



given by the superposition of two independent standing waves in the x and y directions, 

i>(k w ,ky) = A(k x )A(k y ) sin(k x ) sm(k y ), and (22) 

E(k x ,k y ) = -2t(cos(k x ) + cos(k y )). (23) 

In Eqs. ( f22"|) and ( p3[) the pseudomomenta /c x . and k y of a finite L = L x x L y lattice take 
the discrete values k x ( y ) = £ with the integers n.^) = 1, 2, ^(j/)- The normalization 
constants are A(k x ^) = 2/^j2L x{y) + 1 [1 - sin[(2L x( j / ) + l)k x ( y )]/[(2L x(y) + 1) sin k x{y) ]Y l/2 . 
We refer to the quantum numbers k x (y) as the pseudomomenta because the translational 
invariance is broken in a finite lattice, and the momentum is not a good quantum number. 
At = 0.5 0o, the problem reduces to solving the algebraic equation A m (E) = 0, where 
A m for m > 2 obeys a recursion relation A m = [E/t — 2 cos(k y )] A m _! — A m __ 2 with A x = 
(E/t) 2 - 4cos 2 (^) and A = E/t - 2cos(k y ). We note that at = 0.5 O , Eq. (|l|) has 



the symmetries of a bipartite lattice, and therefore the eigenspectrum is highly degenerate 
at this value of the flux. 

At = 0, we can characterize the spectrum using elementary group theory. A finite 
square L = \J~L x \J~L lattice has the symmetries of the C± v spatial point group. It is 
invariant under the following symmetry operations: (a) a rotation of the whole lattice by 
2tt; (b) a rotation by ir; (c) a rotation by ±7r/2; (d) a reflection around the vertical or 
horizontal axes of symmetry; and finally (e) a reflection about the two main diagonals of the 
lattice. These symmetry operations form 5 classes, and therefore the eigenstates of Eq. (|T6|) 
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belong to 5 irreducible representations. The degeneracies are immediately deduced from the 
character tablell of the group C± v given in Table 1. The spectrum is either singly-degenerate, 
belonging to the one-dimensional Al, A2, Bl, or B2 representation, or doubly-degenerate 
belonging to the two-dimensional E representation. In the last row of Table 1 we show the 
characters of the reducible representation of the group C± v of a 3 x 3 lattice in the occupation 
basis. Using the character table of the irreducible representations and the characters of this 
reducible representation, we deduce that for the 3x3 system, group theory predicts 2 pairs of 
doubly-degenerate eigenvalues and 5 singly-degenerate eigenvalues. Similarly, the number of 
degenerate states is deduced for a y/~L x \fL lattice, and the results are given in Table 2. An 
additional electron-hole symmetry in the \J~L x vZ problem mixes the \fL singly-degenerate 
states with total pseudomomenta k x + k y = it at zero energy. This can also be seen from Eq. 
(p3|). To summarize, for a finite square \/L x \[L lattice at = 0, we find g = 1 and g = 2 
degenerate states in the spectrum as well as g = \[~L degenerate states in the middle of the 
spectrum (see Figures 2(a) and 4 at (f) = for the examples of these spectra). As <fi — > 
the magnetic field effects can be calculated in perturbation theory, where the perturbing 
Hamiltonian is given by the total persistent current operator 

H x = -^27T— J2 J2 x [ c \x,y+i) c (x,v) ~ h - c -\- (24) 

<P0 X= l y=l 

At a small flux, the persistent current carried by the singly-degenerate states is linear in 
the flux in the lowest nonvanishing (2nd) order perturbation theory. The persistent current 
carried by the doubly-degenerate states is easily determined by a degenerate first order 
perturbation theory, where the 1st order correction to the energy is E l = ±|(7|i?i|/3) with 
7 and (3 being the two doubly-degenerate states with the same energy at 4> = 0. 

We summarize our results for the limiting values of flux for the 3x3 system in Table 3. 
For generic values of the flux, the spectrum of the 3x3 system is given by the solution of the 
eigenvalue problem of a 9 x 9 matrix. A symmetry operation which commutes with the finite 
two-dimensional square lattice Hamiltonian at <p ^ will reduce the problem. In a particular 
gauge such a symmetry operation should leave invariant the Peirels phase factors in Eq. 
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(PjH. In both the symmetric (A = ^(y, —x, 0)) and the Landau (A = B(0, x, 0)) gauge, upon 
transversing the distance between the two neighboring points (xi,yi) and (x 2 , 2/2) on a square 
lattice, the electron's wavefunction gains a phase proportional to f(2, 1) = 2/2^2 —Xiy\. The 
object /(2, 1) is invariant under a v operations, which are reflections about the main diagonals 
of the lattice, i.e. x — > y, y — > x, and x — > — y, y — > —x. It is easy to verify that this 
symmetry is indeed observed in the numerically calculated current distributions in Fig. 5. In 
these gauges the a v operations exhaust all the spatial symmetries, but we can not rigorously 
prove that this result is gauge invariant. We can, however, comment on the degeneracies 
(or level crossings) which occur in the spectrum in Fig. 2(a) at = | n0 o , n — 1, 2, 3, 4. 
The limiting cases are given in Table 3. The degeneracy at = 0.25 0o at zero energy 
follows from the particle-hole symmetry. The values of flux at which the degeneracies occur, 
= I 0o an d = § 0o ; are the values when the first and second semiclassical orbits are 
commensurate with the lattice. 

In a one- dimensional ring of length L enclosing a magnetic flux 0, the symmetry operator 
which commutes with the Hamiltonian is the magnetic translation operator along the ring, 
and the one-dimensional problem is easily diagonalized for all values of the flux with a 
spectrum given byS 



E n = — 2tcos 



2vr / > 
— I — +n 

'0 j 



(25) 



L 

where n = 1,...,L. At a finite flux, the eigenstates of the ring Hamiltonian are also the 
eigenstates of the momentum, which are obtained from the momentum at = by the 
same shift x(^) ^ or an eigenstates. (As we shall show later, this result also holds in an 
interacting one-dimensional Hubbard ring.) The energy level crossings between any two 
different n and n' states occur at satisfying n — n' = 2-^ + 2 (integer). The resultant 
energy spectrum is periodic in (with a period 0o) through the ring. The discontinuities in 
the persistent current occur at = (0 = 0.5) for even (odd) total number of electrons in 
the ring.0 A typical spectrum of a ID tight-binding ring is shown in Fig. 7 for the L = 9 
site ring. 
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It is of interest to compare the relative magnitudes of the ID and 2D non-interacting 
persistent currents. In the ID ring with N electrons the total current I(N, 0) is of the same 
order as the persistent current carried by the occupied individual states. This happens due 
to the cancelation of the persistent currents carried by the states with opposite momenta in 
the ground state distribution. In one period, I(N, 0) in the ID ring is given by 



2ix \ sin (irN/L) , , . 

I(N odd ) = -J„ sin — f- } '' -0.5 < f - < 0.5, and 26 
\L 0o / sm(7r/L) O 



I(N e 



'• 27r f AT , n <f>\ • /2vr 0\ sin((A^ + 1)tt/L) 
sin — N/2 + f -sm' 1 



< -f < 1, (27) 

00 



L \ 0o / \ L 0o / sin(7r/L) 

where Jo = cvf/ L with vf = 2t/h. Similar formulas were derived in Ref. 0. In a 2D system, 
the magnitude of the total current can decrease due to the cancelation between opposite sign 
persistent current contributions coming from bulk and edge states. In our earlier workjlll we 
showed that the typical 2D persistent current (J 2 ) 1 / 2 scales with the size of the boundary of 
the system. For completeness, we reproduce in Fig. 8 our results for the calculated system- 
size dependence of the typical current in the half-filled ID and 2D systems with the same 
flux through the areas of each system. In the case of a constant flux density (i.e., the same 
flux piercing a unit cell for different systems), the magnetization density (or the persistent 
current) saturates for large system sizes of 2D lattices for all electron fillings as shown in 
Fie. 9. 



VI. THE INTERACTING SPECTRA 
A. The minimal Hubbard model 

1. Periodic S z oscillations 

A puzzling and interesting feature of the results plotted in Figs. 1 and 2 is the periodic 
oscillations of the z-component S z of the total spin in the ground state of a finite quantum 
dot array. The existence of spin flips in the ground state by itself is not surprising. Any level 
crossing between states with different spins leads necessarily to a spin flip. The presence of 
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another flux value within a magnetic period that leads to a reverse flip can not be assumed 
a priori. This implies that the lowest energy states with different S z can be degenerate for a 
range of the magnetic flux or at discrete flux values. We show in Fig. 10 the magnetic flux 
dependence of the ground state energy of the 3x3 array in the minimal Hubbard model 
(Eq. (^)) for values of N and M for which we find the S z oscillations in Fig. 2(b). To obtain 
the results in Fig. 10, we set = = 0, and also leave out the terms not contributing to 

L 

the S z oscillations and consider an interaction of the local form Vi nt = yUpitPu in Eq- ©• 

i=i 

The ground state energy with fixed N and M as a function of the flux appears to consist of 
segments of parabolas with their centers shifted along both the field and energy axes. For 
N = 2, 3, 4, 5, 6 in Fig. 10, the energy parabolas belonging to different M values overlap 
for a range of flux, leading to periodic S z oscillations in Fig. 2(b) for a non-zero Zeeman 
splitting. This result is not restricted to 2D systems. We find it to be valid also in ID rings. 
We show in Fig. 11 the flux dependence of the lowest energy with different N and M for 
which we find S z oscillations in Fig. 2(c) in the Hubbard ring with L = 9 sites. In the rest 
of this subsection, we solve the Bethe ansatz equations to study the S z degeneracy of the 
ground state of a ID Hubbard ring and then generalize our results to two dimensions. 

The energy and the canonical momentum of an L-site Hubbard ring with N electrons 
enclosing a magnetic flux are given by the Bethe ansatz solution 

N 

E = -2 ^cosfy, (28) 

3=1 

N T 9<tt /A 1 9^r I 



^ = E 

3=1 



2vr^ 
3 L<P 



2tt 

T 



(29) 



(For the sake of brevity, we refer to the canonical momentum defined in Eq. (^9j) as the 
momentum from now on.) The energy and the momentum of the interacting system in 
Eqs. (|28|) and ( p9|) are given by expressions similar to those for a non-interacting system, 
involving a summation over N pseudomomenta kj. The pseudomomenta kj describe the 
charge degrees of freedom that are coupled to the spin degrees of freedom with associated 
spin quantum numbers X a in an interacting system. The set of numbers kj and X a , referred 
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to as charge and spin rapidities, is found!! by solving a set of coupled Bethe ansatz equations: 

sin kj — \/3 



M 



Lkj = 2-kI.j + 2tt-J- - 2 tan_1 



N 

2 tan" 1 

3=1 



A a — sin kj 

ujlt 



M 



27T J a + ^ 2 tan 

0=1 



U/At 

A a — 
U/2t 



(30) 
(31) 



The quantum numbers Ij (J a ) are integers if M is even (N-M is odd) and half-odd integers 
if M is odd (N-M is even), and J a is restricted^ to a range \J a \ < (N — M + l)/2. The 
ground state energy is obtained by taking the consecutive sets of integers Ij and J a , and the 
lowest excitations are obtained by having holes in the ground state distribution of quantum 
numbers. 

The problem of the persistent current in a ID Hubbard ring was studied earlier in Ref. 



36| . These authors found that the system accommodates magnetic flux by creating a magnon 
excitation with a hole in the ground state X a distribution. In Ref. ^ it was assumed that 
the spin excitations caused by the magnetic flux in the system remain spin waves for all 
values of U/t. This allowed one0 to conclude that for infinite U / 1, the ground state energy 
has iV cusps in a magnetic period. For large but finite U/t relevant to quantum dot arrays 
under investigation, the assumption that a magnetic field leads to a spin wave excitation 
in a ID ring works well. In Appendix A, we reanalyze the ground state quantum number 
distribution of a Hubbard ring with particular values of N and M that were singled out in 



Ref. |36|. In addition we analyze the ground state distributions of M — 1 states to understand 
the origin of the periodic S z oscillations. 

To summarize our analysis of Appendix A, we find that a ID Hubbard system of TV 
interacting electrons behaves as a single particle on a ring in a magnetic flux: the ground 
state corresponds to a sequence of states with the consecutive values of the total momentum 
which are defined by Eq. (|29|) . The resulting ground state energy has iV intersecting 



parabolic segments per flux period for large U/t in agreement with Ref. [36|. The interacting 
system changes its total momentum as a function of the flux by creating a magnon excitation. 
We find from the numerical solution of the Bethe ansatz Eqs. (p£f ) - ([H|) that the dynamics 
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of this excitation (i.e. its location in the spectrum of spin rapidities) is determined by the 
dynamics of the total momentum that changes sequentially in multiples of 2tt/L from its 
minimum to its maximum value within one magnetic period. This last finding was not 
emphasized in the earlier analysis^ of the persistent current in a Hubbard ring. This result 
follows from Eqs. Q2"8|) - Q3Tj) in the U/N — > oo limit. In the limit of the infinite interaction, 
the pseudomomenta are given by 

2tt 



kj L 



(32) 



A sequential set of the total momentum states which minimizes charge rapidities in Eq. fl32|) 
and, therefore, the ground state energy can be chosen by considering both, either spin or 
charge, excitations in the ground state quantum number distribution. The energy cost to 
create a charge excitation in the infinite interaction limit in a ID Hubbard ring is associated 
with the motion of noninteracting spinless fermions and its magnitude is of the order of t. 
The energy cost to create a spin excitation is determined by t 2 /U, leading to spin excitations 
being the most energy-efficient way in the ID Hubbard ring to accommodate the enclosed 
magnetic flux. 

The total spin S may not change in a new total momentum ground state. Therefore, 
we find that all cusps in the ground state energy within a magnetic period are associated 
with changes in the orbital quantum numbers, but not necessarily with changes of the total 
spin. We find the ground states with nonzero total spin S to be the generic situation with 
respect to the magnitude of the interaction U/t, the electron filling, and the flux. Such 
ground states maintain their (2S + 1) degeneracy with respect to the different values of S z 
over the range of the magnetic flux. This extended degeneracy is maintained until the total 
spin and momentum change in either of the different S z ground states. In the presence of 
the extended total spin degeneracies a small nonzero Zeeman term (u z ) in the Hamiltonian 
in Eq. ([I]) leads to the periodic oscillations of the z-component S z of the total spin S as 
shown in Figs. 1 and 2. The reason is the following. For small uo z ^ the ground state 
is the M = N/2 + 1 state, whereas for to z = the M = N/2 and M = N/2 + 1 states 
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are degenerate. Beyond this region, the M = N/2 state is the lowest energy state, and 
the system changes its S z until it completes a period and the ground state again becomes 
M = N/2 + 1 state. In the 2D 3x3 array, similar extended total spin degeneracies with 
respect to the M, M — 1, M — 2 states lead to the S z oscillations in Figs. 1 (a), (b) and 
2(b). 

2. Spin polarization transition 

The critical magnetic field B c characterizing the spin polarization transition in the min- 
imal Hubbard model is determined to be gfiBB c ~ 7.5^ and « 11^ at half-filling in the 
L = 9 ring and the 3x3 array, respectively (see insets to Figure 1). In an infinite system, 
the critical field B c corresponds to a field driven ferromagnetic transition and was calcu- 
lated by several authors for a ID Hubbard ring neglecting the orbital contributionll. For 
large U/t (strong interaction) the critical magnetic field in the ring depends on the filling 
n = N/L through the relation gnsB c = ^j-n[l — sin(27m)/(27rn)]. In the thermodynamic 
limit, the critical field at half-filling is g\isB c = Qualitatively, the high value (~ llfj) 
of the critical field in the finite 2D system is not surprising, since we expect the quantum 
fluctuations to have a larger disordering effect in 2D systems, raising the critical field. More 
quantitatively, this result can be understood as follows. The critical field determines the gap 
for the triplet excitations in the ferromagnetic phase. At half-filling, the critical field is then 
determined by gfx B B c = E (N = L, M = N - 1) - E (N = L,M = N) with E being the 
lowest eigenvalue of Eq. (JJ) for particular values of iV and M. The energy of the spin po- 
larized state is trivial, E (N, M = N) = —NgfiBB c /2. In a strongly interacting M = N — 1 
case, a spin down electron will 'attract' a spin up hole. The energetically favorable config- 
uration on a lattice corresponds to a single site being occupied by the spin down electron 
and the spin up hole, and all the other sites being singly occupied by the ferromagnetically 
aligned spin up electrons. All the electrons contribute the trivial Zeeman term to the total 
energy, but the site with the spin down electron and its spin up neighbors also contribute 



24 



the exchange energy E exc h. The total energy of the M = N — 1 spin up and one spin down 
electrons is E (N = L, M = N - 1) = -(N - l)g/j lB B c /2 + gfi B B c /2 + E exch . Therefore, the 
critical field is determined by the exchange energy, g/i#.B c = E exc h. The exchange energy is 
given by E exc h = ^j- f(d) with f(d) being the number of spin up neighbors of the spin down 
electron. In a ID ring f(d) = 2, and in large 2D systems on the average f(d) = 4 in the 
nearest neighbor tight-binding square lattice. Note that our numerical result for B c in finite 
ID rings (B c « 7.5t 2 /U) is reasonably close to the thermodynamic result (B c = 8t 2 /U). 



3. Metal-insulator transition and spin ordering 

The charge stiffness, D c , which characterizes [D c = (7^ 0) characterizes the Mott 
insulator (metal)] the degree of localization in a Mott-Hubbard system is defined for a ID 
ring of L sites to 

Note that, by definition, the charge stiffness is the negative flux derivative of the persistent 
current: D c ~ In Fig. 13 we compare the noninteracting and interacting persistent 
currents calculated in the minimal Hubbard model for the 3x3 array. The interacting 
persistent current is suppressed by orders of magnitude at half-filling, n = N/(2L) = 0.5. 
This is a finite size manifestation of the Mott-Hubbard metal- insulator transition as n — > 0.5 
and U/t 7^ 0. For small values of the flux the localization effects due to the magnetic field 
are small, and the magnetic response is not strongly affected by interaction at low filling in 
the metallic phase as can be seen in Figs. 13 (a)-(d). The additional discontinuities of the 
persistent current in Fig. 13 arise from the total spin S or from accidental degeneracies that 
occur with no Zeeman term in the Hamiltonian. In the insulating phase (N = 9 in Fig. 13) 
the persistent current is suppressed by orders of magnitude from its noninteracting value. 
Apparently, the ground state is described by one set of quantum numbers in a magnetic 
period, leading to the absence of discontinuities of the persistent current. The magnetic 
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response is paramagnetic at zero flux. We verify this suppression of the magnitude and 
oscillations as well as the paramagnetic nature of the persistent current at half-filling by 
doing calculations for 2x2, 3x2, 4x2 and 3x3 quantum dot arrays; these results are 
shown in Fig. 14. The rate of the suppression of the persistent current with the increased 
strength of the interaction U/t in the half-filled 3x3 array depends on the orbital magnetic 
flux. We compare this rate in the 3x3 array for the two values of the flux, 0/0o = 0.1 and 
0/0o — 0.4 in Fig. 15. The magnetic field localization at 0/0o = 0.4 enhances the effect of 
the Hubbard U/t on the magnitude of the persistent current, compared to a slower decay of 
the persistent current at 0/0o = 0.1. 

Our use of the terms 'insulating' and 'metallic' to describe the N = L and N < L finite 
Hubbard systems is justified on the basis of the existing workEHll. Using the linear response 
theory it can be shown that the charge stiffness at zero flux in a ID ring is proportional to 
the zero- temperature DC conductivity of the system, with the real part of the conductivity 
given by00 

a(u) = 2tt£) c (0 = 0)6(u) + a reg (u). (34) 

The a reg (ui) term gives the weight of the finite frequency part of the conductivity spectrum, 
and the coefficient D = 2ttD c = iTe 2 (n/m), called the Drude weight, gives the zero frequency 
weight of the conductivity spectrum. In infinite systems, an insulator is characterized by 
D = D c = 0, and a metal has D > oil. However, in finite systems in the insulating phase 
D c is small but nonzero, and is expected to reduce exponentially with increasing system 
size. Stafford and Millis derived an explicit asymptotic form for the scaling function of the 
charge stiffness in a half-filled Hubbard L-site ring (with L even)!!, 

D c (L,0 = 0)U-oo = {-l) L ' 2+1 L^D{U/t)eM-L/aUlt)]. (35) 

In Eq. ( |35"D D{U/t) is a U/t dependent number, and £(U/i) is the localization length, that 
at small U/t behaves asi! 

2t + U/2* + . 

t/-o s A(U,t) v 1 
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In Eq. (^) A(U,t) is the Lieb-Wu Hubbard ga for the charge excitations in the system, 
that is determined by A c (U,t) = E Q (N + 1, M) + E Q {N — 1, M) — 2E (N,M). Stafford 
and Millis also showed^ that in a 1 D Hubbard ring the equal time T = Green's function 
decays exponentially as G(x) ~ exp(— \x\/£). Therefore, the decay of the persistent current 
and that of the single-particle Green's function are shown to be connected in a ID Hubbard 
ring. In the 2D Hubbard model at half-filling, where exact thermodynamic results are not 
available, recent numerical work using Monte-Carlo methods^ demonstrated that the single- 
particle Green's function scales as G(x) ~ exp(— \x\/^2) with a correlation length £2 showing 
a critical behavior different from the ID case. For our purpose the important point is that 
a 2D Hubbard system at half-filling is a Mott insulator, so following the exact ID result of 
Stafford and Millis we make the ansatz that the charge stiffness in a 2D finite system also 
scales in the same way as the 2D Green's function, D c ~ exp —l/^ 2 (U/t) with I ~ being 
the characteristic linear size of the system. This provides some rationale for our assertion 
that N = L results for the finite clusters in Figs. 14 and 15 correspond to the insulating 
phase. 

In the 3x3 array, the persistent current flows along the perimeter of the cluster for 
all N, but the sign of the persistent current in the system can be either diamagnetic or 
paramagnetic. Therefore, a steady paramagnetism that we find in the half-filled small arrays 
can not be understood using a simple analogy with the paramagnetism of the edge states in 
the corresponding noninteracting system!!. The sign of the noninteracting persistent current 
is paramagnetic in the half-filled 2 x 2, 4 x 2, and 3x3 arrays and is diamagnetic in the 
3x2 array. The negative charge stiffness in the interacting 2x2 array is consistent with 
the prediction of Eq. fl35|) for a L = 4 x integer ringll@. We do not know whether this 
orbital paramagnetism is a generic feature of the finite half-filled 2D Hubbard lattices at 
zero temperature, but our finding is consistent with previous exact diagonalization studies 
of the optical conductivity of 4 x 4 periodic 2D Hubbard lattices which also found a negative 
Drude weight at half-fillingS. We have calculated the distribution of the local moments 
((OlS^ilO)) in the half-filled finite lattices with the Zeeman term set equal to zero in the 
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Hamiltonian. The ground state of the 3x3 array belongs to a doublet with the total spin 
S — | and S z = ±|. At U = 0, the uncompensated moment in the system is distributed 
along the diagonals of the lattice. At finite U, the off diagonal dots gain opposite magnetic 
moments, and the system is antiferromagnetically ordered. We verified that the magnitudes 
of the neighboring moments show a tendency to equalize with increasing U/t. Therefore, a 
many-body state in the 3x3 array at half-filling represents a finite-size realization of the 
Mott-Heisenberg insulator^ (antiferromagnetic at a finite U/t) ground state in the system, 
with a Mott-Hubbard gap opening in the charge excitations (as manifested in the strong 
suppression of the persistent current in Figs. 13-15). The ground state for the even L- 
lattices is a singlet (S = 0, S z = 0), with all the local moments being very small. In the ID 
L = 9-site ring (S = \, S z = ±|), we do not find any ordering of the local moments. We 
conclude that in small clean 2D systems the onsite interactions strongly suppress (possibly 
exponentially at iV = L) the magnetic response as N — > L, and for N < L the onsite 
interactions have a much weaker destructive effect. 



B. The extended Hubbard model 

We consider first onsite and nearest neighbor off-site interactions in Eq. (fj|). We set 
Vu = U, and = V if (ij) are nearest neighbor dots (Vij = otherwise). Eq. (H) then has 
the form of an extended Hubbard mo 

The Hubbard U by itself leads to a spin density wave (SDW) commensurate with the 
lattice periodicity at half-filling in a 2D lattice. This SDW state has uniform onsite charge 
density, {p^ = 1. It is easy to see from Eq. ([!]) that at half filling the off-site interaction V 
prefers double occupancy on the alternating sites, e.g. (pi a dd) = 2, and (pi even ) = 0. The 
minimum energy configuration for a large nonzero V is the charge density wave (CDW) state 
commensurate with the lattice. The competition between U and V leads to a SDW-CDW 
transition in both ID and 2D half-filled systenJiHl as V is increased. The mean field phase 
diagram is easily obtained by considering the strong coupling limits of U and V in Eq. ([[]). 
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The energy of the SDW state in the half-filled system is E$dw = + ^N nn with N nn 
being the number of all pairs of nearest neighbors on a lattice. The corresponding energy 
of the CDW state is Eqdw — The phase transition line at half-filling in the extended 

Hubbard model is consequently given by 

N 

V = V C = U—, (37) 

with V < V c being the SDW state. In the 2D lattices with open boundary conditions, 
N nn = 2L — L x — L y yields for large L, V c ~ and in the ID rings N nn = L with 
V c = U. In the usual Hubbard model with periodic boundary conditions, the corresponding 
equations for the SDW-CDW phase transition are V c — j (2D) and V c = ^ (ID). In the half- 
filled ID finite rings the charge stiffness maximizes at the transition line between SDW and 
CDW states@. Intuitively, this seems plausible because the persistent current is sensitive 
to quantum fluctuations, which are maximized at the phase transition. We find a similar 
maximization at half filling (N = 9) of the magnitude of the persistent current in the 3x3 
2D array as a function of the off-site interaction V/t (Fig. 16). The persistent current is 
strongly suppressed on both sides of the SDW-CDW transition while being enhanced at 
the transition point. The current also changes sign at the transition and is more strongly 
suppressed deep in the CDW phase for large U/t. The critical interaction strength V c = 0.65 
at which the finite size SDW-CDW 'transition' occurs in the 2D array does not depend on 
the orbital magnetic flux (for U/t = 10), which is qualitatively consistent with the mean 
field and strong coupling results. The shape of the peak in Fig. 16, however, depends 
strongly on the orbital magnetic flux - the transition at a flux of 0/0o — 0.4 showing a 
larger enhancement than that at 0/0o — 0.1. 

The SDW-CDW transition in a ID extended Hubbard model is foundS to be a second 
order phase transition for U/t < U c /t = 3, and a first order transition for U/t > U c /t. For 
example, the order parameter of the CDW state p = (l/N))X—l) l (f>i) is nonzero in both 
SDW and CDW phases for U > t/ c @S. A recent mean field calculation of the SDW-CDW 
phase diagram finds the transition to be first order (i.e. discontinuous in the order parameter) 
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for all values of U and In the 3x3 system finite size effects dominate for weak coupling, 
and the persistent current peak position depends on the orbital flux for U/t <3. We have 
calculated the persistent current as a function of V/t for U/t — 3, 5, 7, 10, 20, and 30 in 
the half filled 3x3 array and obtained the finite size SDW-CDW phase diagram, which is 
shown in Fig. 17. The slope of ~ 0.68 of the linear fit to the data in Fig. 17 approximately 
agrees with the value 3/4 predicted within the strong coupling theory. 

With the inclusion of the longer range (beyond nearest neighbors) off-site interactions 
in our calculations, we find that the spin antiferromagnetic order is destroyed in the finite 
3x3 array with the electron density distribution remaining uniform. We set the values 
of Vij using the classical capacitance matrix formalism!!! as explained in Section II, and 
show in Fig. 18 our calculated persistent current at half filling (N = 9) as a function of 
V/t = (U/t)(C/C g ), with U/t = {e 2 /C g )/t = 10 being fixed in the 3 x 3 array. An important 
conclusion from Fig. 18 is that strong long range interactions could enhance the magnitude 
of the persistent current back to its non-interacting half filled value, thereby effectively 
negating the Mott-Hubbard localization effect. Whether this is a general result (i.e. valid 
even in the thermodynamic limit) or purely a finite size effect is not clear at this stage. The 
intermediate values of the long-range interaction do not significantly modify the persistent 
current for N < L. This can be seen in Fig. 19 where we compare the {U/t = 10, V = 0) 
and (U/t = 10, V/t — 5) persistent current for various filling fractions in the 3x3 array. 
Our results can be understood on the basis of the dependence of the inverse capacitance 
matrix elements on the relative ratio of C/C g . For simplicity, let us discuss this dependence 
in a double-dot array. In the limit C/C g <C 1, the interdot and intradot interaction matrix 
elements are V±2 = U{C/C g ) — > and Vji = U. For C/C g ^> 1, these matrix elements 
are V 12 = U/2 and V u = (C7/2)(l + C g /C) -> U/2. The large ratio of C/C g in the latter 
case effectively describes a single composite system. In this limit all elements of the inverse 
capacitance matrix become equal. Thus the electrons become totally uncorrelated, and the 
only effect of interactions is to give an overall charging energy to the system. Our results 
suggest that in a finite system with weakly screened interactions the SDW-CDW transition 

30 



may not occur, but the persistent current could still be enhanced to its noninteracting value. 



VII. DISORDERED QUANTUM DOT ARRAYS 

A. Non-interacting case 

In this section we consider the effect of random disorder on the persistent current in the 
finite 2D quantum dot lattices in the spin polarized regime (see Figs. 2(a), 4 and 7 as exam- 
ples of clean system spectra in this regime). We include disorder in our calculation through 
a spin-independent parameter W that denotes the half-width of a uniform distribution of 
random on-site quantum dot energies centered around A. The random on-site single-particle 
energies are set as 

£i = A + Set, 5ei e [-W/2, W/2] . (38) 

The introduction of disorder could, in principle, lead to Anderson localization in the 
system, with all electronic states being exponentially localized in the presence of strong 
disordei@@. In general, disorder (W) and interaction (U, V) compete in determining the 
electronic properties of the resultant (W, U, V ^ 0) Mott-Hubbard-Anderson model. 
Disorder introduces two related length scales in the problem: a mean free path l(W) and 
a localization length £,{W). (In the absence of any interaction, in two dimensions for a 
weak disorder the localization length £(W) depends exponentially^ on the mean free path 
as £ = I explnkp 1/2} with kp being the Fermi wavevector whereas in one dimension £ ~ I.) 
Scaling and perturbative arguments predict that for weak disorder (and in the absence of 
interaction) the conductance of a 2D lattice of size (L) 2 depends logarithmically on the 
mean free path whereas for strong disorder it falls off exponentially with the system size as 
cr(L) ~ exp(— We find two similar weak and strong disorder regimes in the behavior 
of the persistent current as a function of the disorder strength W, indicating a connection 
between the persistent current and the DC conductance of the finite 2D system. In one 
dimension Anderson localization occurs (i.e. er(L) ~ exp(— in the presence of any 

31 



finite disorder^, and the persistent current amplitude is exponentially suppressed^ for all 
W. In Fig. 20 we show the calculated log- log plot at half-filling of the rms current (J 2 ) 1 / 2 , 
averaged over 100 disorder realizations for each value of W, as a function of the disorder 
strength W for various array sizes (3x3, 4x4, 5x5, 6x6). In plotting these results, we 
have factored out a scale factor, n c / L = (L 1 / 2 — l) 2 /L, so that the results for various system 
sizes fall on top of each other, showing approximate scaling with system size and disorder. 
(The scale factor n c = (L 1 / 2 — l) 2 is the number of unit cells or plaquettes in each square 
array of size L.) The two dashed straight lines in Fig. 20 give the best fits to weak and 
strong disorder scaled currents, leading to the following empirical results for the effect of 
disorder on the persistent current: (I 2 )^ = (L/n c ) g(W), with the scaling function g(W) 
being given by 

[ w-i ,7 = (6.4 ±2.8) x 10~ 2 , W< 1.55 7rf 
g(W)~\ (39) 

I W~P, 13 = 1.84 ±0.49, W> 1.55 7rt. 

We note that the empirical scaling defined by Eq. ^ (and shown in Fig. 20) is consistent 
with well-known noninteracting scaling localization result of there being logarithmic 2D 
localization at weak disorder. 



B. Interacting case 

Random disorder effects on the persistent current are subtle in the presence of onsite 
interactions in the 3x3 array: they depend both on the filling in the array and the relative 
ratio U/W. Away from half-filling, the random disorder suppresses the magnitude of the 
current and smoothens its discontinuities arising from level crossings of energies with differ- 
ent orbital quantum numbers. The S z oscillations and the corresponding discontinuities in 
the current in the metallic phase persist in the presence of intermediate ( W ~ U) disorder. 
Hubbard interaction effects are not particularly important at high disorder, which is what 
is seen in Fig. 20 (asterisks) for the typical interacting current even at half filling. The in- 
termediate disorder (W < U) produces an antilocalization effect at half-filling by enhancing 
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the persistent current from its finite ^/-suppressed value as can be seen in Fig. 21 (a) and 
(b). We also observe that the disordered persistent current is enhanced in the presence of the 
long range interaction from its [/-suppressed value as can be seen in Fig. 21(c). Finally we 
note that the presence of disorder does not reverse the roles of the weak and strong orbital 
flux: the relative behaviors of persistent current at 0/0o = 0.1 and 0/0o — 0.4 in Fig. 21 
can be compared with those for clean systems plotted in Figs. 15 and 18. 



VIII. CONCLUSIONS 

We calculate the persistent current and electron addition spectrum in coherent two- 
dimensional semiconductor quantum dot arrays and one-dimensional quantum dot rings 
pierced by an external magnetic flux, using the exact diagonalization and the Bethe ansatz 
techniques within an extended Mott-Hubbard Hamiltonian. We find that the magnetization 
density of a finite multidot array is periodic in the magnetic flux. In weak fields, we find 
flux periodic oscillations in the S z component of the total spin S. We have included in our 
model the effects of both intradot/interdot Coulomb interactions and random disorder. We 
find that the persistent current is suppressed in the antiferromagnetic Mott-insulating phase. 
The finite-size realization of the spin density wave - charge density wave ordering transition 
has been found to maximize the 2D array persistent current at half-filling at the critical value 
of the nearest neighbor interaction, a behavior qualitatively similar to the charge stiffness of 
ID rings0. We obtain the phase diagram for the SDW-CDW transition. We demonstrate 
that the electrostatic long-range interdot interactions enhance the magnitude of the Mott- 
insulating phase persistent current to its non-interacting system value. We find that the 
noninteracting persistent current as a function of the random disorder strength exhibits a 
behavior similar to that of the conductivity^ of the system: it is strongly suppressed only at 
large disorder strengths. The Anderson-Mott transition^ has been found to have a subtle 
effect on the persistent current in the 3x3 array at half filling: the intradot Coulomb 
interaction less than or comparable to the disorder strength W (U ~ W) increases the 
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disordered system persistent current. 

We believe that the 2D persistent current physics in the 3x3 array discussed in this 
paper has already been indirectly observed in transport measurements in a 3 x 3 array of 
strongly coupled quantum dots fabricated by means of square grid gate structures on top 
of a GaAs/AlGaAs heterostructureS. In low magnetic fields {B < 1 T) at T rs 40 mK the 
three conductance minima were measured in Ref. | at B=0.18 T, B=0.48 T, and B=0.80 
T, and superimposed on these minima small oscillations periodic in B with period of 12 
niT, 21 niT, and 24 mT were seen! in the measured magnetoconductance of the 3x3 array. 
The large conductance minima were interpreted to be due to the classical localized orbits 
in the array, with the cyclotron radius of each orbit being equal to R c = m*vp/eB, with 
R c = 330 nm, 125 nm, and 75 nm in the three minima, respectively. The measured period 
of small oscillations was attributed to a reduction of a single classical orbit Aharonov-Bohm 
period due to the coupling between all possible classical orbits with a given R c in the lattice. 
We suggest here an alternative explanation of this experiment. First, the observed large 
conductance minima correspond to the first three Landau orbits being commensurate with 
the size of a single dot. The semiclassical cyclotron radii, R n = J (Ti/eB){2n + 1), are 60.5 
nm (n=0), 64 nm (n=l), and 64 nm (n=2) for the three conductance minima, respectively. 
Therefore, we estimate the radius of each dot as 63 nm (which is reasonably consistent with 
the other dimensions of the array). These minima should be observed for temperatures 
for which hu c < 2kT yielding T « 1.8 if, which is consistent with T = 2 K where the 
large conductance minima start to be affected by the temperature!. We agree with the 
experimentalists that the shorter period oscillations are due to the coupling of the electron 
orbits in the lattice. But we think that the contributing orbits originate from the quantum 
mechanically required gauge invariance in the experimental phase coherent system. These 
orbits arise due to electron hopping between the localized states of each dot. We estimate 
that the magnetic flux corresponding to one flux quantum 0q piercing a unit cell of the lattice 
is ~ 46 mT (with a lattice constant of a ~ 300 nm). For the density of the 2DEG and 
the estimated size of the dot, the number of electrons in each dot is ~ 5. Considering one 
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mobile electron at the Fermi energy in each dot, the experimental 3x3 array approximately 
corresponds to our 3x3 array at half filling (in the situation of the weakly screened long- 
range Coulomb interactions). The persistent current in our case is found to oscillate with a 
period of 0/0o = 0.25 and 0.5 (see Figs. 13(h) and 19(h)), that corresponds to AB = 11.5 
mT and AB = 23 mT, which are consistent with the experimentally measured values. 

Finally, we emphasize three essential limitations of our work which may (generally) re- 
strict its quantitative applicability to realistic quantum dot arrays (we do believe that our 
results describe well the qualitative aspects of coherent and collective physics in semicon- 
ductor quantum dot arrays): (1) we ignore completely the complicated (and interesting) 
details of single dot electronic structure, approximating the energy level in each dot by two 
spin-split energy levels described by the simple Fock-Darwin-Zeeman model; (2) we use the 
Mott-Hubbard model in treating interaction and correlation effects (both intradot and inter- 
dot) - such a simple parameterization of Coulomb correlations may not apply quantitatively 
in real quantum dot arrays; (3) we have considered only small (no more than 9 dots) finite 
arrays, being limited entirely by the exponentially growing Hilbert space size in our extended 
Hubbard model. While we see no particular hope of going beyond these approximations in 
the near future, we should emphasize that these limitations do not substantially restrict 
the qualitative applicability of our results to the coherent collective physics of quantum dot 
arrays. For example, currently fabricated coherent quantum dot arrays typically contain 
only two to four quantum dots, and therefore our finite size calculations are, in fact, per- 
fectly appropriate. Also, the basic physics of Coulomb blockade and quantum fluctuations 
are entirely captured in our extended Mott-Hubbard model, and the complications of single 
dot energetics (beyond that of two-level Fock-Darwin-Zeeman physics) can usually be dis- 
tinguished from the collective physics of interest to us. Finally, some of our non-essential 
approximations, for example, restricting to zero temperature and using a constant interdot 
hopping amplitude t, can easily be relaxed in future calculations if such a need arises. We 
emphasize here that these non-essential approximations do not affect our qualitative conclu- 
sions in any significant manner. We believe that the extended finite Mott-Hubbard model as 
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used by us is the appropriate minimal model which should form the basis of discussing the 
collective and coherent physics in mesoscopic semiconductor quantum dot arrays. We feel 
that the direct observation of equilibrium persistent current and spin oscillations in coherent 
quantum dot square lattices and rings will shed light on the interplay among the single dot 
physics, coherence, disorder, long- and short-range Coulomb interaction effects on quantum 
phase transitions that we discuss in this paper. 
This work was supported by the U.S.-O.N.R. 

APPENDIX A: GROUND STATE ENERGY OF A HUBBARD RING 
ENCLOSING A MAGNETIC FLUX 

Here we analyze the ground state quantum number distribution of a Hubbard ring that 
we obtain by solving numerically the Bethe ansatz Eqs. ( p8| ) - (|3l| ) for all N < L, M = N/2 
and M = N/2 - 1. 

a. N = An + 2, M = 2n + 1, and M = 2n. For N = An + 2 and M = 2n + 1 the ground 
state distribution at a flux <p > is given byS: 

{J . } = _(^_ 1 ,_(^ _(^ +p _ 1 ,_(^ +p+1 ,...,(^. (A2) 



In Eq. ( |A"2l ) the position p of a hole in the distribution of the spin rapidities varies from p = 
to p = M within one half of a magnetic period, and p = corresponds to the ground state 
distribution at <ft = 0. The ground state distribution for negative values of the flux is obtained 
by shifting all J a to the right by one unit. The momentum P of the state is given by the 
momentum q of the spin wave excitation P = q = — ^ p. Thus, the total momentum is zero 
at zero flux and then takes the values of all consecutive multiples of — ^ within one half of 
magnetic period, with the momentum being equal to -^P = —M = — y at (p/cpo = 0.5. The 
distribution in Eq. is valid for N < L, for N = L the distribution of both the charge 



and spin quantum numbers remains symmetric at a nonzero flux. Since -k-P m ax — m 
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the limit of large U/t the ground state energy has N cusps at 0/0o = ±(p + 1) and the 
persistent current is diamagnetic at 0/0o — ^Pjn- It was proved by Stafford and Millisii 
that the N = 4n + 2 ground state of the electron system under periodic boundary conditions 
is a spin singlet (total spin S — 0) at — 0, and it can be either a spin triplet (total spin 
S = 1) or a singlet at = 0.5 0o (a flux value where the antiperiodic boundary conditions 
are realized). In the L = 9-site Hubbard ring we find that at = 0, the ground state is a 
spin singlet, and at = 0.5 0o it is a spin triplet with S = 1. At = 0.5 0o, the ground 
state, therefore, has a (2S + 1) spin degeneracy. 

The energy of the N = An + 2 and M = 2n (S z = 1) state is higher than the ground 
state energy of the singlet state at = (as an example, see Fig. 11(d) for iV = 6, M = 3, 
and M = 2). The energy of the S z = 1 state at > is minimized by choosing the ground 
state distribution as 

N N N 
{Ij} = ~ y, — g + 1, y - 1, and (A3) 

M = JM^, _ p , _ „ + 2 , H^i) + , (A4) 

In Eq. ( |A4p the position p of a hole in the distribution of the spin rapidities varies from 
p = M to p = within one half of a magnetic period, and p = M corresponds to a 
ground state distribution at <ft = 0. The charge degrees of freedom contain a nonzero 
momentum = — y in the system, and the spin rapidities are positioned in the Fermi sea 
to minimize the total momentum at = 0. The vacant hole in this system is also within 
the Fermi bounds. The total momentum of the state described by Eqs. ( |A3|) and (|A4]) is 
t^tP = — (y — p). Thus, at = and 4>/4>q = 0.5, the minimum and maximum momenta 
are = —1 and = — y, respectively. Within our numerical accuracy we do not 
find that the excitation can create a zero total momentum at <fi = in the range of the 
interaction U/t values from to 200. Therefore, at = 0, there is a cusp corresponding to 
the crossing of equal and opposite minimum momenta, and the ground state energy has a 
total of N — 1 cusps per magnetic period. The spin quantum numbers of the N = L ground 
state remain symmetric in a magnetic period. The S z = and the S z — 1 states have the 
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same total momenta at the last segments of the energy curves that includes = 0.5 0o- 
This is a situation for the N = 6, M = 3 and N = Q, M = 2 states in the L = 9-site 
ring at = 0.5 0o in Fig. 11(d). The total spin of the M = 2 system does not change 
in a magnetic period and is equal to 1. The degenerate M = 3 and M = 2 states belong 
to a spin triplet in the segment of the energy curve centered at = 0.5 O . We conclude 
that the magnetic flux changes the total momentum of a state sequentially by one unit from 
one parabola segment to another (for smaller values of U/t, particular values of the total 
momenta are missing in the ground stated). Each energy parabola segment has a fixed 
value of the total spin, which may or may not change with a change of the total momenta. 
The extended degenerate nonzero total spin regions are not specific to the ground state of 
the N = An + 2 systems. They occur whenever the magnetic field drives the systems with 
different S z through a series of the consecutive total momentum states to the flux region 
where the ground state has a nonzero total spin 5 and the different S z states thus have the 
same total momenta. 

b. N = An, M = 2n, and M = 2n - 1. For N = An (N < L) and M = 2n the ground 
state distribution at a flux > is the same as for N = An + 2 and M = 2n case that is 
given in Eqs. ( |X3] ) and (|A4|). In the current situation, a spin wave excitation is above the 



Fermi sea, and the minimum momentum is = — — M) = at (f> = 0. The ground 
state energy has iV cusps in a period,il and the persistent current is diamagnetic at the 
(>o = ±pj^. The ground state distribution for N = An and M = 2n — 1 is given in Eqs. 



lT|) and (|A2|) . For large U/t, an excitation above the negative Fermi sea of spin rapidities 
can be created, so the momentum at 0/0o — 0.5 is = — y. The S z = and S z = 1 
states can have the same minimum momenta P = at <fi = 0. In Fig. 11(c) for N = A 
electrons in the ring, the energy parabolas of M = 2 and M = 1 belonging to the spin 
triplet (5=1) states are coincident in the P = magnetic flux region centered at <fi = 0. 
For iV = 8 in Fig. 11(e) the levels of the P = momentum triplet state (5 = 1, S z — 1) 
and j^P = ±4 singlet state (5 = 0, S z = 0) cross in the vicinity of = and are non 
degenerate at = 0. 
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c. N = An + I, M = 2n, and M = 2n - 1. For N = An + 1 (N < L) and M = 2n 
electrons in a Hubbard ring, both sets of Ij and J a are integer numbers. To accommodate 
the variation of the total momenta, two holes are present in the ground state distribution 
at a nonzero positive flux: 

(N-l) (N-l) (N-l) , lk , 

{Ij} = -^"2"^. -^-^ + !. and ( A5 ) 

(M — 2) M M M 

{«U = - 2 »-»-y +P2~2,-y +p 2 ,...,0,...,p 1 -l, Pl + l,...,y. (A6) 



In Eq. ( |A6|) the position pi of a hole in a distribution of spin rapidities varies from 
Pi = at <p — to pi = y with p 2 being fixed at p 2 = 0. While the hole p\ transverses 
from to y, the momentum varies from to — 4^. A further increase of the momentum 
is accomplished through the motion of the second hole p 2 from to y - 1 with p\ being 
fixed at 4f- The total momentum becomes P = —^-M at 4>/4>q = 0.5. Therefore, in the 
ground state distribution there is always one hole which moves. We show the energies for 
consecutive values of the total momentum in a L = 15 ring for n = 2 and L = 9 ring for 
n — 1 for U/t — 200 in Fig. 12 (a) and (b), respectively. The ground state energy curves 
in Fig. 12 (a) and (b) consist of the N consecutive momentum states within a magnetic 
period. The earlier analysis of the persistent current in a Hubbard ring did not emphasize 
the dynamics of the total momentum, and Yu and Fowler mistakenly concluded that their 
Eq. (3.13) describes the ground state distribution for any nonzero flux, and that at the one 
half flux quantum all the J a 's have to be consecutive integers^. 

The charge and spin quantum numbers are half-odd integers for N = An + 1 and M = 

2n — 1 electrons in a ring. The ground state distribution at <fi > for n > 1 is given by 

N N N 
{I J } = -j,-j + l,-,j-l,and (A7) 

r _ M M M 13 1 Ad . . 

Ua} = -y,.-,-y +p 2 -l,-y +P2 + 1, Pi,-, y + 1. (A8) 



In Eq. ( |A8|) the hole p\ moves from to (M + 3)/2 with p 2 being fixed at p 2 = 0, and then 
the hole p 2 moves from to (M — 3)/2 while p\ = (M + 3)/2. The motion of the holes 
in the distribution of the spin rapidities leads to a cancelation of the —N/2 momentum 
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concentrated in the charge degrees of freedom at zero flux and the consequent integral 
decrease of the total momentum. In the L = 9 ring with N = 5 electrons, for the value of 
the interaction U/t = 10, the states M = 2 and M = 1 have different values of the total 
spin S = \ and S — § and the total momenta j^P = — 1 and -^P = —2, respectively. These 
states are nondegenerate, and there are no oscillations of the total spin S z in the ground 
state of the L = 9, N = 5 Hubbard ring. 

d. N = 4n — l,M = 2n — l, and M = 2n — 2. The lowest energy state for N = 4n — l, 
M = 2n — 1 is obtained by choosing the distribution of quantum numbers that is given in Eqs. 
(|A7|)- (|A8|) . At zero flux, the hole pi starts moving from J a = |, and it moves (M + l)/2 — 1 
consecutive steps to the right. In the remaining flux region the hole P2 transverses to the 
right in (M — l)/2 — 1 consecutive steps. The distribution in the M = 2n — 2 case is 
given in Eqs. flAlp - flAlf) . In this case the hole P2 moves until the total momentum becomes 
— (A^ + l)/2 at 0/0o = 0.5. At = for = 3 electrons in a L = 9-sites ring in Fig. 11(b), 
the two lowest M = 1 and M = energy states with momentum P = and the total spin 
S — I are degenerate in the = region. 
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FIGURES 

FIG. 1. (a) [(c)] A chemical potential /ijy as a function of flux <p in units of 4>q = h/e through 
an elementary cell [a ring in (c)] of a 3 x 3 quantum dot lattice [9— site ring in (c)] with N{= 1 — 9) 
excess quasiparticles in the minimal Hubbard model approximation; (b) [(d)] A corresponding 
z-component S z of the ground state total spin as a function of 4>. The magnetic flux is rescaled by 
1/32 to show the entire dynamics of the spectrum on one flux scale as explained in the text. Insets: 
The critical magnetic flux (f> c of the spin polarization transition as a function of the electron filling 
n = N/ (2L) in the arrays. 

FIG. 2. (a) A chemical potential spectrum plotted versus flux through a unit cell of a 
non-interacting 3x3 array with £\ = A, N = M in Eq. (1); (b) [(c)] A chemical potential 
as a function of flux (top panel), and the corresponding S z component of the ground state total 
spin (bottom panel) of a 3 x 3 array [9— site ring in (c)] with = A, e-\ = A(l + <5/A) (S/A = 0.03) 
in Eq. (1). (a), (b) and (c) can be directly compared with high and low field regimes in Figs. 1(a) 
and (c), respectively. 

FIG. 3. A spectrum of Eq. (16) plotted for rational fractions of flux <p/4>Q = p/q (q = 120, p 
is incremented from 1 to 60) through a unit cell of an infinite tight-binding lattice. The first two 
Landau bands (n = 0, 1) are marked on the plot. 

FIG. 4. A spectrum of a finite 15 x 15 tight-binding lattice versus 4>/4>o through a unit cell. 
The different regions that are marked (a)-(l) on the plot are discussed in the text. 

FIG. 5. The persistent current and charge density distribution of the single-particle states 
marked in Fig. 4. The darkness of circles, and darkness, or thickness of connecting lines are 
proportional to the magnitude of charge and current, respectively. Arrows on the lines show a 
direction of the persistent current. The maximum current and charge in (a) to (1) in units of (et/h) 
and e, respectively, are 0.028, 0.009; 0.114, 0.018; 0.232, 0.023; 0.250, 0.102; 0.393, 0.039; 0.150, 
0.030; 0.162, 0.062; 0.172, 0.019; 0.126, 0.011; 0.161, 0.027; 0.073, 0.015; 0.099, 0.017. 
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FIG. 6. The total persistent current I calculated using Eq. (15) for (a) N = 26 electrons on a 
15 x 15 lattice((e) in Figs. 4 and 5), and (b) N = 73 ((h) in Figs. 4 and 5). 

FIG. 7. The energy spectrum of a 9-site tight-binding ring shown for one period of flux through 
the area of a ring; this spectrum can be compared with the high field regime in Fig. 1(c). 

FIG. 8. {I 2 ) 1 / 2 at half filling versus system sizes (L) (L = L 1 / 2 x L 1 / 2 with L 1 / 2 varying from 
2 to 10): triangles (2D array), squares (ID ring). (The average is taken over the range of the same 
total flux through the systems, which varies from —(f>o/2 to 4>o/2.) The dashed lines are the best 
linear fits to the dependence {I 2 ) 1 / 2 ~ L a with a = 0.46 (1.1) in a 2D (ID) systems. 

FIG. 9. A typical current (I 2 ) 1 / 2 , I at 4> = 0.25 <fio, and the maximum current I max are shown 
from top to bottom, respectively, as a function of N electrons in L 1 / 2 x L 1 / 2 systems (with L 1 / 2 
varying from 2 to 10) with the same flux density for systems of different sizes. 

FIG. 10. The energy of the 3x3 quantum dot array calculated using Eq. (1) in a minimal 
Hubbard model approximation (^p- = U = Wt, Vij = 0, e-f = e± = 0) as a function of the flux 
through the unit cell. The energy is plotted for all values of N and M which are found in the 
ground state of the array in Fig. 2(b). In (a)-(f), a number of electrons iV in the array is as 
indicated in the legend; M = N/2, M = N/2 — 1, M = N/2 — 2 energy states are shown with the 
solid, dotted, and dashed lines, respectively. The states with different M that are degenerate over 
the range of flux have the total spin in (a) S = 1; (b) S = §; (c) 5 = 2, 1, 0; (d) S = §; (e) S = 1. 
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FIG. 11. The energy of an 9— site Hubbard ring as a function of the flux through the ring 
calculated using the Bethe ansatz Eqs. (28)-(31). All parameters and plotting conventions are 
as given in the caption in Fig. 10. A sequence of states in one magnetic period, for the flux 
-0.5 < 0/0o < 0.5 is given by (a) ±P = 1,0,-1 (M = 1) and ±P = 1,-1 (M = 0); (b) 
±P = 1,0,-1 (M =1) and ^P = 0(M = 0); (c) ^P = 2,0,-2 (M = 2) and Ap= 1,0,-1 
(M = 1); (d) ±P = 3, 0, -3 (M = 3), and = 3, 1, -1, -3 (M = 2); (e) ^P = 4, -4 (M = 4) 
and ^P = (M = 3). The states with different M that are degenerate over the range of flux have 
the total spin in (a) S = 1; (b) S = §; (c) S = 1; (d) S = 1. 

FIG. 12. The energy of the ID Hubbard rings versus the flux is shown for N = 4n+ 1, M = 2n 
case for L = 15, A = 9, M = 4, 17 = 200 t (a) and L = 9, A = 5, M = 2, [7 = 200 t (b) for the 
consecutive values of the total momentum in the rings. Each energy curve is labeled in the plot by 
the value of its total momentum expressed in units of 7^; in (b) it is also labeled by the value of 
the total spin S. 

FIG. 13. The persistent current in the 3x3 array of quantum dots as a function of flux through 
a unit cell for U/t = (solid line) and U/t = 10 (dashed-dotted line) for A = 2 - 9 in (a)-(h), 
respectively. The array is modeled in the minimal Hubbard approximation, with other parameters 
set the same as in Fig. 2(b). The corresponding z-component S z of the total spin in a ground state 
is shown in the bottom part of the plots, with an arbitrary offset for the illustrative purposes. 

FIG. 14. The persistent current at half-filling in the 2x2, 3x2, 4x2 clusters of quantum 
dots with S = 0, S z = ground state and in the 3x3 lattice with S = \, S z = \ ground state 
as a function of flux through the unit cell. A key to each curve is indicated in the legend. All the 
parameters set the same as in Fig. 2(b). 

FIG. 15. The persistent current as a function of the Hubbard U ~/t in the 3x3 array with 
iV = L = 9 at (f)/4>Q = 0.1 (solid) and 4>/4>o = 0.4 (dotted). The same parameters used as in Fig. 
2(b). 



47 



FIG. 16. The persistent current as a function of the nearest neighbor off-site interaction V/t in 
the 3 x 3 array with N = L = 9 and U/t = 10 at <p/(p = 0.1 (solid) and <j)/(p = 0.4 (dotted). The 
same parameters used as in Fig. 2(b). 

FIG. 17. A phase diagram of the finite size SDW-CDW transition in the 3x3 array with 
N = L = 9 is shown for the onsite and off-site nearest neighbor interaction U/t and V/t at 
4>/4>o = 0.1 (squares) and (f)/4>o = 0.4 (triangles). The same parameters used as in Fig. 2(b). 

FIG. 18. The persistent current as a function of the long range off-site interaction 
V/t = (U/t){C/C g ) in the 3 x 3 array with N = L = 9 and U/t = 10 at 0/>o = 0.1 (solid) 
and 4>/4>o = 0.4 (dotted). The same parameters used as in Fig. 2(b). 

FIG. 19. The persistent current in the 3x3 array of quantum dots as a function of flux through 
a unit cell for U/t = 10 and long range interaction V/t = 5 (solid line) and U/t = 10, V/t = 
(dashed-dotted line) for N = 2 — 9 in (a)-(h), respectively. Other parameters are set the same as in 
Fig. 2(b). The corresponding S^-component of the ground state total spin is shown in the bottom 
part of the plots, with an arbitrary offset for the illustrative purposes. 

FIG. 20. A log-log plot of (/ 2 ) 1/2 versus W averaged over 100 disorder realizations for each 
value of W is shown (U = V = 0) at half-filling in 3 x 3 (crosses), 4x4 (triangles), 5x5 (solid 
line), 6x6 (pluses). Asterisks show the interacting (U/t = 10, V = 0) results at half-filling for the 
3x3 system. 

FIG. 21. The persistent current as a function of (a) and (b) the Hubbard U/t, (c) the long 
range off-site interaction V/t = (U/t)(C/C g ) (with U/t = 10) in a disordered 3x3 array with 
N = L = 9 and at (f)/4>o = 0.1 (solid) and 4>/4>o = 0.4 (dotted). The persistent current is averaged 
over 10 disorder realizations. The disorder strength W is indicated in the plots. The same other 
parameters used as in Fig. 2(b). 
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TABLE I. Character table of C± v taken from Ref. [27]. The last row gives characters of C4, 
on a 3 x 3 lattice in the occupation basis. 
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TABLE II. A number m(i?) of eigenvalues belonging to an irreducible representation R of the 
C4D group on a 

L l/2 x L l/2 lattice _ 
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TABLE III. A summary of the analysis of the 3x3 tight-binding lattice spectrum for the 
limiting values of the flux. 
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Figure 1. 
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Figure 3. 

R. Kotlyar, C. A. Stafford, and S. Das Sarma, Addition spectrum, persistent current 




Figure 4. 
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Figure 5. Continued 

R. Kotlyar, C. A. Stafford, and S. Das Sarma, Addition spectrum, persistent current .... 




Figure 6. 

R. Kotlyar, C. A. Stafford, and S. Das Sarma, Addition spectrum, persistent current 




Figure 7. 

R. Kotlyar, C. A. Stafford, and S. Das Sarma, Addition spectrum, persistent current 




Figure 8. 

R. Kotlyar, C. A. Stafford, and S. Das Sarma, Addition spectrum, persistent current 




Figure 9. 

R. Kotlyar, C. A. Stafford, and S. Das Sarma, Addition spectrum, persistent current 



(a) N = 2 (b) N=3 




Figure 10. 
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